WEEK 11 NOTES
#Bike Share Demand Forecasting with Panel Data & Temporal Lags
The dreaded empty bikeshare station:
Operational question:
“At 6:00 AM, which stations will run out of bikes by 8:00 AM?”
Your job: Build a system that predicts demand in space AND time
Previous weeks (Weeks 5-7):
This week:
Definition: Data that follows the same units over multiple time periods
Cross-sectional data:
Panel data:
| Station ID | Date-Hour | Trip Count | Temperature | Day of Week |
|---|---|---|---|---|
| 1 | 2018-05-01 08:00 | 12 | 65°F | Tuesday |
| 1 | 2018-05-01 09:00 | 15 | 67°F | Tuesday |
| 1 | 2018-05-01 10:00 | 8 | 69°F | Tuesday |
| 2 | 2018-05-01 08:00 | 5 | 65°F | Tuesday |
| 2 | 2018-05-01 09:00 | 7 | 67°F | Tuesday |
Each row = station-hour observation with features and outcome
Raw trip data:
Trip 1: Started at 8:05:23 AM
Trip 2: Started at 8:07:41 AM
Trip 3: Started at 8:15:12 AM
Trip 4: Started at 8:23:08 AMHourly binning:
# All trips between 8:00-8:59 AM → "08:00" bin
dat <- dat %>%
mutate(interval60 = floor_date(ymd_hms(start_time), unit = "hour"))Result:
Now we can count: “Station A had 15 trips in the 8:00 AM hour”
Finer temporal resolution:
dat <- dat %>%
mutate(interval15 = floor_date(ymd_hms(start_time), unit = "15 mins"))Trade-offs:
Today: We’ll use hourly bins for simplicity
dat <- dat %>%
mutate(
interval60 = floor_date(ymd_hms(start_time), unit = "hour"),
week = week(interval60), # Week of year (1-52)
dotw = wday(interval60, label=TRUE), # Day of week (Mon, Tue, ...)
hour = hour(interval60) # Hour of day (0-23)
)These become predictors:
hour %in% c(7,8,9, 17,18,19)dotw %in% c("Sat", "Sun")Core idea: Past demand predicts future demand
Spatial features (Week 6):
Temporal features (This week):
Intuition: If there were 15 trips at 8 AM, there will probably be ~15 trips at 9 AM
study.panel <- study.panel %>%
arrange(from_station_id, interval60) %>% # Sort by station, then time
group_by(from_station_id) %>%
mutate(
lag1Hour = lag(Trip_Count, 1), # Previous hour
lag2Hours = lag(Trip_Count, 2), # 2 hours ago
lag3Hours = lag(Trip_Count, 3), # 3 hours ago
lag12Hours = lag(Trip_Count, 12), # 12 hours ago
lag1day = lag(Trip_Count, 24) # Yesterday same time
) %>%
ungroup()Important: Lags are calculated WITHIN each station
Station A on Monday:
| Time | Trip Count | lag1Hour | lag12Hours | lag1day |
|---|---|---|---|---|
| Mon 07:00 | 5 | NA | NA | NA |
| Mon 08:00 | 12 | 5 | NA | NA |
| Mon 09:00 | 15 | 12 | NA | NA |
| Mon 19:00 | 8 | 10 | 5 | NA |
| Tue 08:00 | 14 | 6 | 10 | 12 |
Interpretation: At Tue 08:00, demand was 14. Same time yesterday (Mon 08:00) was 12.
Different lags capture different patterns:
Model will learn which lags are most predictive for each station/time combination
Remember Week 6?
# Spatial features we created:
boston.sf <- boston.sf %>%
mutate(
crimes_500ft = ..., # Count crimes nearby
crime_nn3 = ..., # Average crime at 3 nearest neighbors
dist_downtown = ... # Distance to downtown
)This week is the same concept:
# Temporal features we're creating:
study.panel <- study.panel %>%
mutate(
lag1Hour = ..., # Demand at nearby TIME
lag1day = ..., # Demand at similar TIME
)Problem: Not every station has trips every hour
# Data we have (sparse):
Station_A, May_1_08:00, 12_trips ✓
Station_A, May_1_09:00, 0_trips ✗ (missing row!)
Station_A, May_1_10:00, 8_trips ✓But we NEED:
# Complete panel (every station-hour combination):
Station_A, May_1_08:00, 12_trips
Station_A, May_1_09:00, 0_trips ← Must exist with 0!
Station_A, May_1_10:00, 8_tripsWhy? Lag calculations break if rows are missing
Step 1: Calculate all possible combinations
# How many unique stations?
length(unique(dat_census$from_station_id)) # e.g., 600 stations
# How many unique hours?
length(unique(dat_census$interval60)) # e.g., 744 hours (31 days)
# Total combinations needed:
600 stations × 744 hours = 446,400 rowsReality check: Do we have 446,400 rows? Probably not!
# Create every possible station-hour combination
study.panel <- expand.grid(
interval60 = unique(dat_census$interval60),
from_station_id = unique(dat_census$from_station_id)
)
# Join to actual trip counts
study.panel <- study.panel %>%
left_join(
dat_census %>%
group_by(interval60, from_station_id) %>%
summarize(Trip_Count = n()),
by = c("interval60", "from_station_id")
) %>%
mutate(Trip_Count = replace_na(Trip_Count, 0)) # Fill missing with 0Now every station-hour exists, even if Trip_Count = 0
Each station has fixed characteristics:
# Station location, demographics from census
station_data <- dat_census %>%
group_by(from_station_id) %>%
summarize(
from_latitude = first(from_latitude),
from_longitude = first(from_longitude),
Med_Inc = first(Med_Inc),
Percent_White = first(Percent_White),
# ... other demographics
)
# Join to panel
study.panel <- study.panel %>%
left_join(station_data, by = "from_station_id")Result: Every row has station location + demographics
Some features change over time:
# Weather changes hourly
weather <- weather_data %>%
select(interval60, Temperature, Precipitation)
study.panel <- study.panel %>%
left_join(weather, by = "interval60")
# Create time features
study.panel <- study.panel %>%
mutate(
week = week(interval60),
dotw = wday(interval60, label = TRUE),
hour = hour(interval60),
weekend = ifelse(dotw %in% c("Sat", "Sun"), 1, 0)
)What we have now:
Ready to model! But first… temporal validation
Remember when we talked about spatial lag models?
Spatial Lag Model (SAR):
\[Y_i = \rho \sum_j w_{ij}Y_j + \beta X_i + \varepsilon\]
Why NOT for prediction:
Our Approach Today:
\[Y_{it} = \alpha_i + \beta_1 X_{it} + \beta_2 Y_{i,t-1} + \varepsilon_{it}\]
Station fixed effects + temporal lags:
Key difference: We use PAST outcomes as features, not NEIGHBOR outcomes
You CANNOT train on the future to predict the past!
WRONG approach:
# Train on Weeks 3-4
train <- data %>%
filter(week >= 19)
# Test on Weeks 1-2
test <- data %>%
filter(week < 19)This is predicting the past using the future!
** CORRECT approach:**
# Train on Weeks 1-2
train <- data %>%
filter(week < 19)
# Test on Weeks 3-4
test <- data %>%
filter(week >= 19)This is predicting the future using the past!
Real-world scenario:
It’s May 15, 2018. You need to forecast demand for May 16-31.
# Split by time (week of year)
train <- study.panel %>%
filter(week < 19) # Weeks 1-2 of May (early period)
test <- study.panel %>%
filter(week >= 19) # Weeks 3-4 of May (later period)
# Fit models on training data only
model <- lm(Trip_Count ~ lag1Hour + lag1day + Temperature + weekend,
data = train)
# Evaluate on test data
predictions <- predict(model, newdata = test)This mirrors operational deployment: Predict tomorrow using yesterday’s patterns
Week 7: Spatial cross-validation
This week: Temporal validation
Both are about out-of-sample generalization!
We’ll build 5 models, adding complexity:
Goal: See which features improve prediction accuracy
model1 <- lm(
Trip_Count ~ hour + dotw + Temperature + Precipitation,
data = train
)Predictors:
hour: Hour of day (0-23)dotw: Day of week (Mon, Tue, …, Sun)Temperature: Hourly temperaturePrecipitation: Rainfall amountCaptures: Daily and weekly cycles, weather effects
model2 <- lm(
Trip_Count ~ hour + dotw + Temperature + Precipitation +
lag1Hour + lag3Hours + lag1day,
data = train
)New predictors:
lag1Hour: Demand 1 hour agolag3Hours: Demand 3 hours agolag1day: Demand 24 hours agoHypothesis: Past demand predicts future demand
model3 <- lm(
Trip_Count ~ hour + dotw + Temperature + Precipitation +
lag1Hour + lag3Hours + lag1day +
Med_Inc + Percent_Taking_Public_Trans + Percent_White,
data = train
)New predictors:
Captures: Neighborhood effects on demand
model4 <- lm(
Trip_Count ~ hour + dotw + Temperature + Precipitation +
lag1Hour + lag3Hours + lag1day +
Med_Inc + Percent_Taking_Public_Trans + Percent_White +
as.factor(from_station_id), # Station dummies
data = train
)Station fixed effects:
Why? Some stations are just busier than others!
model5 <- lm(
Trip_Count ~ hour + dotw + Temperature + Precipitation +
lag1Hour + lag3Hours + lag1day +
Med_Inc + Percent_Taking_Public_Trans + Percent_White +
as.factor(from_station_id) +
holiday + holiday_lag1 + holiday_lag2, # Holiday indicators
data = train
)Memorial Day weekend (May 28-29, 2018):
holiday: Current hour is during holidayholiday_lag1: 1 hour into holidayholiday_lag2: 2 hours into holidayCaptures: Disruption to normal patterns during holidays
Mean Absolute Error:
\[MAE = \frac{1}{n}\sum_{i=1}^{n}|y_i - \hat{y}_i|\]
Interpretation:
Example: MAE = 5 means predictions are typically off by ±5 trips
| Model | MAE (Test Set) | What Improved? |
|---|---|---|
| 1. Baseline | 8.2 trips | — |
| 2. + Lags | 6.5 trips | Temporal persistence |
| 3. + Demographics | 6.1 trips | Neighborhood effects |
| 4. + Fixed Effects | 5.3 trips | Station baselines |
| 5. + Holidays | 5.1 trips | Disruption events |
Biggest improvement: Adding temporal lags (8.2 → 6.5)
Not all errors are equal!
Question for operations: When do prediction errors cause bikes to run out?
Map errors by station:
When are we most wrong?
Scatterplot by time of day:
Diagonal red line: Perfect predictions
Blue line: Our actual performance
What you might see:
Critical question: Are errors related to demographics? (Equity concern!)
# Join errors back to demographic data
station_errors <- station_errors %>%
left_join(
dat_census %>%
distinct(from_station_id, Med_Inc,
Percent_Taking_Public_Trans, Percent_White),
by = "from_station_id"
)
# Plot relationships
station_errors %>%
pivot_longer(cols = c(Med_Inc, Percent_Taking_Public_Trans,
Percent_White),
names_to = "variable", values_to = "value") %>%
ggplot(aes(x = value, y = MAE)) +
geom_point(alpha = 0.4) +
geom_smooth(method = "lm", se = FALSE) +
facet_wrap(~variable, scales = "free_x")For a bike rebalancing system:
Who benefits from accurate predictions?
Higher accuracy in:
Lower accuracy in:
Critical analysis: Are we reinforcing unequal service quality?
Questions to consider:
What could reduce errors?
weekend * hour)Common policy applications:
Anytime you have: Same units observed repeatedly over time
Your task: Build space-time models for Chicago bike share
Working through the code:
expand.grid() (15 min)Discussion questions:
Key functions you’ll use:
# 1. Binning time
floor_date(datetime, unit = "hour")
# 2. Creating complete panel
expand.grid(station = stations, time = times)
# 3. Temporal lags
group_by(station_id) %>% mutate(lag1 = lag(outcome, 1))
# 4. Fixed effects
lm(y ~ x + as.factor(station_id))
# 5. Evaluation
MAE = mean(abs(observed - predicted))Open the lab file: week11_bike_share_lab.Rmd
Work in groups of 2-3
I’ll circulate to help with:
Take your time - this is complex material!